Methods and systems for determining gas permeability of a subsurface formation

ABSTRACT

Methods and systems for determining permeability, as a function of pore pressure, and porosity of a subsurface formation. The method includes positioning a sample in a sample assembly comprising of a gas and a pressure gauge, inside a pressure vessel comprising gas or liquid and a pressure gauge, measuring a first gas pressure, p i , of the sample inside the pressure vessel, applying a second gas pressure, p o , to the pressure vessel, the second gas pressure being greater than the first gas pressure, measuring a third gas pressure, p, at time, t, at location, x, from the inlet of sample inside the pressure vessel, determining a total gas mass per unit volume of the subsurface formation, m, and determining the permeability, k, of the subsurface formation as a function of pore pressure based at least in part on the first gas pressure, the second pressure, the third gas pressure, and the gas density, with a single test run.

CROSS-REFERENCE TO RELATED PATENT APPLICATIONS

This application claims priority from U.S. patent application Ser. No.15/378,305, filed Dec. 14, 2016, which claims priority from ProvisionalApplication No. 62/267,091, filed Dec. 14, 2015, both titled “METHODSAND SYSTEMS FOR DETERMINING PORE PRESSURE-DEPENDENT GAS PERMEABILITY OFA SUBSURFACE FORMATION,” the entire contents of which is incorporatedherein by reference.

TECHNICAL FIELD

Example embodiments relate to methods and systems for determiningpermeability and porosity of a subsurface rock formation using pore gaspressure.

BACKGROUND

Unlike conventional reservoirs, pores in shale formations are extremelysmall, typically on the order of nanometers. In these nano pores, anon-negligible portion of gas molecules collides more often with thepore wall than with other molecules, and thus so-called “slip flow” andKnudsen diffusion occur. Previous studies on gas flow in shale matrixfound that the gas permeability in shale is a function of the pore gaspressure because the slip flow and Knudsen diffusion effect becomessignificant when the pore gas pressure is relatively low.

Shale gas permeability as a function of pore gas pressure, resultingfrom “slip flow” and diffusion processes, is critical for characterizingand modeling gas flow in a shale gas reservoir. However, this importantpore gas pressure-dependency is hardly considered in practice because ofthe lack of a practical and efficient technique that can be usedroutinely for determining the pressure-dependent shale gas permeability.

Pressure dependence has a significant impact on predicted gas-productionrate. There are currently two approaches to measure the pressuredependence of gas permeability in the laboratory. The first one is tosimply perform a number of pulse-decay permeability tests underdifferent gas pressures. Then, these tests will provide gas permeabilityvalues for a number of gas pressures. Initially, the system is inequilibrium with a given gas pressure. A small pressure pulse is thenintroduced into the upstream gas reservoir, such that the pulse does nothave a significant disturbance to the gas pressure in the system. Thepressures at the two gas reservoirs are monitored as a function of time.The pressure evolution results are fitted using analytical solutions,with permeability being a fitting parameter. However, it generally takesa relatively long time to equilibrate the test system from one testpressure to the next one.

The other approach to determine the pressure dependence is to firstdevelop a formulation of gas permeability as a function of gas pressureand then estimate values for parameters in the formulation bynumerically matching the relevant test results under different gaspressure conditions. Test results are generally different frompulse-decay tests in which the pressure pulse is not limited to a smallone because numerical model is flexible enough to incorporate the pulsedisturbance to the system. However, non-uniqueness of parameterestimation is always a problem for inverse modeling. Also the accuracyof estimated results from this approach is ultimately determined by thatof the used formulation of gas permeability as a function of gaspressure that is not fully established yet.

SUMMARY

Example embodiments disclosed provide a new method to measurerelationship between shale gas permeability and pore gas pressure. Thedevelopment is based on a new analytical solution to one-dimensional gasflow under certain boundary and initial conditions. The advantages ofthe disclosed approach over the currently available ones include that itdirectly measures the relationship using a single test run and withoutany presumption regarding the form of parametric relationship betweengas permeability and pressure. In addition, the current approach allowsfor estimating both shale permeability and porosity at the same timefrom the related measurements.

One example embodiment is a transient flow method for determiningpermeability of a subsurface formation. The method includes extracting asample of the subsurface formation, positioning the sample in a pressurevessel comprising a type of natural or man-made gas or liquid and apressure gauge, measuring a first pore gas pressure, p_(i), or theinitial pore gas pressure of the sample inside the pressure vessel,applying a second pore gas pressure, p_(o), to the inlet of the samplewithin pressure vessel, the second pore gas pressure being greater thanthe first pore gas pressure (p_(o)>p_(i)), measuring a third pore gaspressure, p, at location x as a function of time t along the sample inthe pressure vessel, determining a total gas mass per unit volume of thesubsurface formation, m, and, determining the permeability function,k(p), (hereinafter referred to as permeability, and k) of the subsurfaceformation from gas transport parameter D(p), based at least in part onthe first pore gas pressure, the second pressure, the third pore gaspressure as a function of time, and the gas density, with a single testrun. The relationship between D(p) and permeability is given in Equation9 below. The method may also include determining the gas transportparameter of the subsurface formation, D(p), using a first formula:

${D(p)} = {- \frac{\int_{p_{i}}^{p}{\frac{\lambda}{2}\frac{dm}{dp}{dp}}}{\frac{dp}{d\;\lambda}}}$

where p_(i) is the first pore gas pressure inside the sample before thesecond pore gas pressure p_(o) is applied, p is the third pore gaspressure at location x at time t, m is the gas density or total gas massper unit volume of the subsurface formation, and λ is an independentvariable calculated using the formula xt^(−1/2). Then permeability maybe determined from D(p) using Equation 9.

The method may also include determining the total gas mass per unitvolume of the subsurface formation, m, using a second formula:m=ϕρ+(1−ϕ)ρ_(a)

where ϕ is porosity of the subsurface formation, ρ is gas density of thefree gas, and ρ_(a) is adsorbed gas mass per unit volume of thesubsurface formation.

The method may also include determining the porosity ϕ of the subsurfaceformation using a third formula:

$\phi = \frac{B - {A{\int_{p_{i}}^{p_{0}}{\lambda\;\frac{d\;\rho_{a}}{dp}{dp}}}}}{A{\int_{p_{i}}^{p_{0}}{\lambda\;\frac{d\left( {\rho - \rho_{a}} \right)}{dp}{dp}}}}$

where A is a cross-sectional area of the sample, and B is a slope of acurve of the cumulative gas flow into the sample at x=0 (sample inlet)versus t^(1/2).

The method may also include determining the slope of the curve, B, usinga fourth formula:

$B = {A{\int_{p_{i}}^{p_{0}}{\lambda\;\frac{dm}{dp}{dp}}}}$

Another example embodiment is a non-transitory computer-readable mediumhaving computer executable instructions that cause a computer to performthe operations of reading a measurement of a first pore gas pressure,p_(i), of a sample inside a pressure vessel comprising a sample assemblyof a subsurface formation, a type of natural or man-made gas or liquid,and a pressure gauge, reading a measurement of a second pore gaspressure, p_(o), applied to the inlet of the sample, the second pore gaspressure being greater than the first pore gas pressure, reading ameasurement of a third pore gas pressure, p, as a function of time t, atlocation x from the end of the sample close to the inlet in the pressurevessel (hereinafter referred to “location x”), determining a total gasmass per unit volume of the subsurface formation, m, and determining apermeability of the subsurface formation, k, based at least in part onthe first pore gas pressure, the third pore gas pressure, and the gasdensity, with a single test run.

The computer executable instructions further cause the computer toperform the operation of determining the transport parameter of thesubsurface formation D(p) using a first formula:

${D(p)} = {- \frac{\int_{p_{i}}^{p}{\frac{\lambda}{2}\frac{dm}{dp}{dp}}}{\frac{dp}{d\;\lambda}}}$

where p_(i) is the first pore gas pressure inside the pressure vesselbefore the second pore gas pressure p_(o) is applied, p is the thirdpore gas pressure at location x at time t, m is the total gas mass perunit volume of the subsurface formation, and λ is an independentvariable calculated using the formula xt^(−1/2). Then permeability maybe determined from D(p) using Equation 9.

The computer executable instructions further cause the computer toperform the operation of determining the total gas mass per unit volumeof the subsurface formation, m, using a second formula:m=ϕρ+(1−ϕ)ρ_(a)

where ϕ is porosity of the subsurface formation, ρ is gas density of thefree gas, and ρ_(a) is adsorbed gas mass per unit volume of thesubsurface formation.

The computer executable instructions further cause the computer toperform the operation of determining the porosity ϕ of the subsurfaceformation using a third formula:

$\phi = \frac{B - {A{\int_{p_{i}}^{p_{0}}{\lambda\;\frac{d\;\rho_{a}}{dp}{dp}}}}}{A{\int_{p_{i}}^{p_{0}}{\lambda\;\frac{d\left( {\rho - \rho_{a}} \right)}{dp}{dp}}}}$

where A is a cross-sectional area of the sample, and B is a slope of acurve of the cumulative gas flow into the sample at x=0 (the sampleinlet) versus t^(1/2).

The computer executable instructions further cause the computer toperform the operation of determining the slope of the curve, B, using afourth formula:

$B = {A{\int_{p_{i}}^{p_{0}}{\lambda\;\frac{dm}{dp}{dp}}}}$

Another example embodiment is an apparatus for determining thepermeability of a subsurface formation. The apparatus includes a firstsleeve having a length, a diameter, a first open end, and a second openend, a first end piece adapted to be inserted into the first open endand a second end piece adapted to be inserted into the second open endof the first sleeve, a first hole formed through the first sleeve andinto a body of the sample, a half sleeve disposed on the first sleeve,the half sleeve having an assembly on a second hole corresponding to thefirst hole on the first sleeve, consisting of a tubing that connects toa pressure gauge, the tubing passing the second hole, inserted throughthe first hole, and into the body of the sample, an anchoring device forsecuring the tubing to the first sleeve and half sleeve, thereby formingan assembly, and a pressure vessel for receiving the assembly, thepressure vessel comprising a fluid and a plurality of pressure taps,wherein at least one of the pressure taps is coupled to the pressuregauge in the aforementioned tubing.

According to one example embodiment, the first sleeve and the halfsleeve comprise at least one of rubber and a polymeric material.According to one example embodiment, an inner diameter of the halfsleeve is slightly smaller than the outer diameter of the first sleeve.According to one example embodiment, a length of the half sleeve isequal to or less than the length of the first sleeve. According to oneexample embodiment, the fluid comprises natural gas, water, or oil.According to one example embodiment, the tubing may also include atemperature gauge. According to one example embodiment, pressure gaugemay be coupled to the pressure tap using a flexible line. According toone example embodiment, the apparatus may further include an inlet pumpconfigured to pump gas from a first gas tank into the pressure vessel,and an outlet pump configured to store gas from the pressure vessel intoa second gas tank.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG. 1 illustrates an example system for determining permeability andporosity of a subsurface formation, according to one example embodimentof the disclosure.

FIGS. 2A-2D illustrate an example apparatus for determining permeabilityand porosity of a subsurface formation, according to one exampleembodiment of the disclosure.

FIG. 3 illustrates an example apparatus for determining permeability andporosity of a subsurface formation, according to one example embodimentof the disclosure.

FIGS. 4A and 4B illustrate example apparatuses for determiningpermeability and porosity of a subsurface formation, according to someexample embodiments of the disclosure.

FIG. 5 shows a flow chart illustrating example operations in a methodfor determining permeability and porosity of a subsurface formation,according to one example embodiment of the disclosure.

FIG. 6 shows a graphic illustrating simulated time values when theboundary effect occurs (actual time) and estimated time at differentlocations along a shale formation sample, according to one exampleembodiment of the disclosure.

FIG. 7 shows a graphic illustrating simulated gas pore gas pressure,real value and with random error of +/−0.2 psi, as a function of time atthe core location about one inch away from the inlet, according to oneexample embodiment of the disclosure.

FIG. 8 shows a graphic illustrating comparison between the truepermeability as a function of pore gas pressure and the permeabilityfrom numerical experiment results, according to one example embodimentof the disclosure.

FIG. 9 shows a programmable computer and various forms of computerreadable media, according to some example embodiments of the disclosure.

DETAILED DESCRIPTION

Turning now to the figures, FIG. 1 illustrates an example system 10 fordetermining permeability function, k(p), (hereinafter referred to aspermeability, and k) and porosity of a subsurface formation, accordingto some example embodiments of the disclosure. System 10 includes asubsurface formation sample assembly 300 (as shown in FIG. 2), insidewhich is a sample 30, such as a shale sample or a limestone sample or asandstone sample, in the form of a cylinder or column (referred to ascolumn hereinafter), that may be extracted from the subsurface fordetermining characteristics of the formation. Sample assembly 300 isintroduced in a pressure vessel 50 that may contain a confining fluid22, such as gas or a water-based fluid or an oil-based fluid. However,example embodiments described herein refer to systems with oil as theconfining fluid. The centerpiece of the measurement system is the sampleassembly 300 connected to the pumps through the pressure lines with allmonitoring and regulating devices; and the pressure vessel 50 providesthe confining pressure to the sample assembly 300 and the interfaces tothe tubing and lines of measurements via various types of couplings.

System 10 includes an inlet pump 16 configured to pump gas from a firstgas tank 12 into the core sample assembly 300, and an outlet pump 18configured to pump gas from a second gas tank 14 into core sampleassembly 300. Both pumps may include precise pressure and flowratecontrol and measurement. Pressure vessel 50 may be equipped with ahydraulic pump 20 that may pump oil 22 into and may include apparatusthat monitor and regulate the pressure within the pressure vessel 50.High accuracy temperature and pressure gauge 34 is connected to rocksample and 32 to the inlet of core sample assembly 300; both of themhaving high accuracy transducers to measure temperature and pressure,respectively in real-time. An inlet 28 to the core sample assembly 300may be diverted at a plurality of points using bypass valves 24 andoutlet pipe 26 in order to regulate the pore gas pressure (e.g., theestablishment of the initial pore pressure) in the sample 30 inside thesample assembly 300 which is placed in the pressure vessel 50; the inlet40 and outlet 42 on the sample assembly 300 (as shown in FIG. 2) beingconnected through couplings on the wall or the end caps of the pressurevessel to the pressure lines (such as 28) that are connected the inletpump 16 and outlet pump 18. The pressure vessel 50 may also be equippedwith additional temperature and pressure gauges such as gauge 34, whichmay be in direct contact with the shale sample 30.

FIGS. 2A-2D illustrate in further detail an example set up for arrangingthe shale sample assembly 300 in the pressure vessel 50. As illustrated,the shale sample 30 is first inserted into a sleeve 52 having a length,a diameter, a first open end 60, and a second open end 70. The shalesample may be enclosed in the sleeve 52 using a first end piece 42adapted to be inserted into the first open end 60 and a second end piece40 adapted to be inserted into the second open end 70 of the sleeve 52.A through hole or a port for pressure measurement 65 is formed throughthe body of the sleeve 52 and into the body of the sample 30 so as toinsert a tubing, such as a tubing 44. A half sleeve 46 may be disposedon the sleeve 52, and the half sleeve may include a second holecorresponding to the first hole 65 on the first sleeve. The tubing 44may include a temperature gauge 34 and a pressure gauge 32 asillustrated in FIG. 1, for example. The tubing 44 is inserted throughthe first hole, the second hole, and into the body of the sample 30. Ananchoring device 48, for example, may be used for securing the tubing 44to the half sleeve 46, thereby forming an assembly. Other fasteningdevices, such as ring clamps, may be used to secure the half sleeve 46and the sleeve 52. FIG. 3 illustrates a cross-sectional view of theapparatus in FIG. 2D where steel tubing 44 is inserted through the halfsleeve 46, sleeve 52, and secured using anchoring devices 48. Accordingto one example embodiment, sleeve 52 and half sleeve 46 may include atleast one of rubber and a polymeric material. According to anotherexample embodiment, an inner diameter of the half sleeve 46 may besmaller than the outer diameter of the sleeve 52. According to anotherexample embodiment, a length of the half sleeve 46 is equal to or lessthan the length of the sleeve 52.

After the sleeve 52 is secured for preventing leakage from the port ofpressure measurement 65, the assembly is disposed in the pressure vessel50, as illustrated in FIG. 4A, for example. The pressure vessel 50 mayinclude a plurality of pressure taps such as 54, each of which to beconnected a pressure measurement location along the rock sample. In thiscase, multiple measurement locations along the rock sample may beinstrumented according to the procedure described in FIGS. 2 and 3. Thetubing may also include a temperature gauge 34. According to one exampleembodiment, pressure gauge 32 connected to 44 may be coupled to one ofthe pressure taps 54 using a flexible line 56.

Analytical Method for Determining Permeability and Porosity of aSubsurface Formation

The following sections provide an example method for determiningpermeability, k, and porosity of a subsurface formation using the system10 illustrated in FIG. 1. The method is based on a new analyticalsolution to one-dimensional gas flow under certain boundary and initialconditions, which will be described in further detail below. Thegoverning mass balance equation for gas flow may be given by Equation 1as follows.

$\begin{matrix}{\frac{\partial m}{\partial t} = {\frac{\partial}{\partial x}\left( {\frac{k\;\rho}{\mu}\frac{\partial p}{\partial x}} \right)}} & (1)\end{matrix}$

where t is time, x is the spatial coordinate (a distance from the inletof the sample along its axis), k is the permeability, μ, ρ, and p aregas viscosity, density and pressure, respectively, (note k, μ, ρ, arefunctions of p) and m is the total gas mass per unit volume of theporous medium or apparent gas density, which may be given by Equation 2as follows.m=ϕρ+(1−ϕ)ρ_(a)  (2)

where ϕ is porosity and ρ_(a) is adsorbed gas mass per unit volume ofsolid phase or the subsurface formation. For conservative gases, thesecond term on the right hand of Equation 2 can be considered to bezero.

In Equation 1, the storage term can be rewritten as

$\begin{matrix}{\frac{\partial m}{\partial t} = {\frac{dm}{dp}\frac{\partial p}{\partial t}}} & (3)\end{matrix}$

The present method may relate to isothermal conditions, and therefore mmay be considered a function of pressure only. The method may alsoinclude relatively high confining stress such that the effect ofmechanical deformation due to pore gas pressure change can be ignored.Accordingly, the contributions of gas density change in pressure tostorage can be given by Equation 4 as follows.

$\begin{matrix}{\frac{dm}{dp} = {{\phi\;\frac{d\;\rho}{dp}} + {\left( {1 - \phi} \right)\frac{d\;\rho_{a}}{dp}}}} & (4)\end{matrix}$

Taking into consideration an infinite long shale sample in the form of acylinder/column with gas flow from the inlet (x=0) and subject to thefollowing boundary and initial conditions:p(x,t)=p _(i) (x≥0,t=0)p(x,t)=p ₀ (x=0,t>0)p(x,t)=p _(i) (x→∞,t>0)  (5)

where p_(i) is the initial pressure inside the measurement system beforethe elevated upstream pressure, p₀, is applied.

Using the transformation

$\begin{matrix}{\lambda = {xt}^{- \frac{1}{2}}} & (6)\end{matrix}$

Equations 5 and 1 can be transformed as follows.

$\begin{matrix}{\begin{matrix}{{p(\lambda)} = p_{i}} & \left( \lambda\rightarrow\infty \right) \\{{p(\lambda)} = p_{0}} & \left( {\lambda = 0} \right)\end{matrix}\mspace{14mu}{and}} & (7) \\{{{- \frac{\lambda}{2}}\frac{dm}{dp}\frac{dp}{d\;\lambda}} = {{\frac{d}{d\;\lambda}\left\lbrack {{D(p)}\frac{dp}{d\;\lambda}} \right\rbrack}\mspace{14mu}{where}}} & (8) \\{{D(p)} = \frac{k\;\rho}{\mu}} & (9)\end{matrix}$

Equation 8 is an ordinary differential equation with a as the onlyindependent variable.

Directly integrating Equation 8 for the interval (λ, ∞) yields

$\begin{matrix}{{D(p)} = {- \frac{\int_{p_{i}}^{p}{\frac{\lambda}{2}\frac{dm}{dp}{dp}}}{\frac{dp}{d\;\lambda}}}} & (10)\end{matrix}$

It indicates that D(p) can be fully determined when p(λ) is known.

Based on the gas mass balance, the cumulative gas flow into the column(at x=0) can be determined using Equation 11 as follows.

$\begin{matrix}{{M(t)} = {{A{\int_{0}^{\infty}{\left( {m - m_{i}} \right){dx}}}} = {{{{A\left( {m - m_{i}} \right)}x}❘_{0}^{\infty}{{- A}{\int_{p_{0}}^{p_{i}}{x\frac{dm}{dp}{dp}}}}} = {A{\int_{p_{i}}^{p_{0}}{x\frac{dm}{dp}{dp}}}}}}} & (11)\end{matrix}$

where A is the cross-sectional area of the shale column. CombiningEquations 11 and 6 gives

$\begin{matrix}{{M(t)} = {{\left( {A{\int_{p_{i}}^{p_{0}}{\lambda\frac{dm}{dp}{dp}}}} \right)t^{\frac{1}{2}}} = {Bt}^{\frac{1}{2}}}} & (12)\end{matrix}$

where B is a slope for the curve of M(t) versus t^(1/2). CombiningEquations 4 and 12 gives

$\begin{matrix}{\phi = \frac{B - {A{\int_{p_{i}}^{p_{0}}{\lambda\frac{d\;\rho_{a}}{dp}{dp}}}}}{A{\int_{p_{i}}^{p_{0}}{\lambda\frac{d\;\left( {\rho - \rho_{a}} \right)}{dp}{dp}}}}} & (13)\end{matrix}$

Equations 10 and 13 may be used for estimating gas permeability andporosity, according to one or more example embodiments of thedisclosure.

As illustrated above, for an infinite long shale column with a uniforminitial pore gas pressure, one can directly estimate the porosity andpermeability as a function of pore gas pressure using Equations 10 and13 from measurement of M(t) and p(λ) obtained under a constant pressureat the column inlet. It should be noted, however, that the gascompressibility and adsorption parameter, which are functions of poregas pressure, in these equations may be independently determined orestimated from other tests or existing literature. The latter may not beinvolved if gas used for a test is not reactive.

The test method is consistent with initial and boundary conditions usedto obtain Equations 10 and 13 and M(t) and p(λ) can be effectively andreliably measured from a test run. Initially, shale column 30 with aconfining stress has a uniform gas pore-pressure p_(i). The confiningstress may be significantly higher than the range of pore gas pressureused in the test such that mechanical deformation due to pore gaspressure variation can be ignored. The column 30 may be about 4″ longsuch that it is long enough to be treated as infinite long for a certainperiod of test time. The upstream gas reservoir of the shale column 30may then be connected to pump 16 with precise pressure and/or flowratecontrol. The upstream pressure of the core column 30 may be maintainedas a constant p₀ by the pump 16. The pressure range between p_(i) and p₀covers the range of practical interest or the range in which thepressure dependence is important. Cumulative gas mass flow rate into thecolumn inlet, M(t), may be monitored. The pore gas pressure may bemeasured as a function of time at a given location of shale column thatcan be set any location except two ends, preferably about 1″ from thecolumn inlet. From the transformation given in Equation 6, p(λ) can beobtained from the pressure measurements. The pressure at outlet of thecolumn 30 is measured to monitor pore gas pressure breakthrough.Pressure breakthrough is considered to occur at the outlet when pressureincreases by about 0.1 psi. It should be noted, however, that afterpressure breakthrough, the boundary effect from the downstream may bepropagated to the measurement point. After that time (t_(c)), the lengthof column 30 cannot be treated as infinite anymore. Thus only pressuredata before that time (t_(c)) can be used to calculate p(λ).

The time t_(c) can be estimated using Equation 14 as follows.

$\begin{matrix}{t_{c} = {t_{b}\left\lbrack {1 + \left( \frac{L_{b}}{L} \right)^{2}} \right\rbrack}} & (14)\end{matrix}$

where t_(b) is the time of the pressure breakthrough at the outlet ofcolumn 30, L is the length of shale column 30, and L_(b) is the distancebetween a pressure measurement location and column outlet. The aboveequation may be obtained by assuming D(p) in Equation 9 to be constant.In this case, the travel distance of diffusion front resulting fromoutlet disturbance may be proportional to the square root of the timesince the pressure breaks through at the outlet.

Example Method for Determining Permeability and Porosity of a SubsurfaceFormation

Turning now to FIG. 5, illustrated is a flow chart showing exampleoperations in a method 500 for determining permeability and porosity ofa subsurface formation, according to one example embodiment of thedisclosure. The method uses only one pressure measurement locationbetween the inlet and the outlet of a shale column. However, this isonly for illustration purposes, and the method 500 may include pressuremeasurement at multiple locations along the length of the shale column.While the theory requires only one location to make the pressuremeasurements as a function of time, pore gas pressures at two or morelocations may be measured for reasons like achieving better resolutionfor p(λ). It should be noted, however, that p(λ) can be constructed withpressure measurements at different locations. At each measurementlocation, a small hole with diameter of 1/16″ or less and with a depthto about the center of the shale column 30 may be drilled such thatshale pore gas pressure can be reliably measured and at the same time, asmall hole may not introduce a considerable disturbance to the gas flowalong the column.

At operation 502, the core sample may be assembled with pressuremeasurement equipment and leakage prevention feature, and the sample isplaced in a pressure vessel, as illustrated in FIG. 1, for example. Atoperation 504, an initial pore gas pressure may be established in theshale column and confining stress may be imposed until the pore gaspressure becomes equilibrium. At operation 506, the pore gas pressure inthe gas reservoir connected to the column inlet may be raised to therequired pressure so that gas may flow through the column. Thedifference (p_(o)−p_(i)) between this pressure and the initial pressurecovers the pressure range of interest. At operation 508, the pressuresmay be measured as a function of time from the shale column and at theoutlet and monitor the cumulative flow rate into the column. At thispoint, the pressure breakthrough time (t_(b)) at the outlet may bedetermined and t_(c) can be calculated using Equation 14 as follows.

$t_{c} = {t_{b}\left\lbrack {1 + \left( \frac{L_{b}}{L} \right)^{2}} \right\rbrack}$

At operation 510, p(λ) may be determined using pressure data at timesmaller than t_(c) and based on Equation 6 as follows.

$\lambda = {xt}^{- \frac{1}{2}}$

At operation 512, D(p) may be determined based on Equation 10 and usingp(λ) obtained from operation 510 as follows.

${D(p)} = \frac{\int_{p_{i}}^{p}{\frac{\lambda}{2}\frac{d\; m}{dp}{dp}}}{\frac{dp}{d\;\lambda}}$

In Equation 10, the derivative dp/dλ can be estimated from dp/dt:

$\begin{matrix}{\frac{dp}{d\;\lambda} = {{- \frac{dp}{dt}}\frac{2t^{3\text{/}2}}{x}}} & (15)\end{matrix}$

It may be more convenient to estimate dp/dt because pressure is directlymeasured at location x as a function of t. The values of dp/dt can beestimated using the finite difference method with time interval of 1second or less.

At operation 512, the permeability, k, can be determined as a functionof pore gas pressure using D(p) and Equation 9 as follows.

${D(p)} = \frac{k\;\rho}{\mu}$

At operation 514, the porosity may be determined with p(λ) obtained fromoperation 508 and Equation 13 as follows.

$\phi = \frac{B - {A{\int_{p_{i}}^{p_{0}}{\lambda\frac{d\;\rho_{a}}{dp}{dp}}}}}{A{\int_{p_{i}}^{p_{0}}{\lambda\frac{d\;\left( {\rho - \rho_{a}} \right)}{dp}{dp}}}}$

The above method is based on an analytical solution for infinitely longcolumn. However, in the experimental data a 4″ long shale column, not afinitely long column, is used, and to make sure that the infinitely longcolumn assumption practically holds, pressure data before time t_(c)given in Equation 14 may be used. FIG. 6 shows a graphic illustratingsimulated time values 600 (based on Equation 1 and related boundary andinitial conditions) when the boundary effect occurs (actual time fromsimulation) 602 and estimated time (from Equation 14) 604 at differentlocations along a shale formation sample, according to one exampleembodiment of the disclosure. The initial pore gas pressure in a samplemay be 100 psi. Pressure at the inlet may be instantaneously raised to1000 psi at t>0. Gas density and viscosity in Equation 1 may be treatedas functions of pore gas pressure. Two columns with lengths of 4″ and12″ are used in the simulations. No pore gas pressure breakthrough isobserved for the long column during the test time periods; thus it canbe treated as an infinitely long column. The simulated pressures atdifferent locations for both columns are compared. The outlet boundaryeffect is considered to occur at a time when the pressure difference forthe two columns at a location is larger than 0.1 psi. FIG. 6 showssimulated time when the boundary effect occurs 602 and the timeestimated 604 from Equation 14 at different locations along the column.Clearly, estimates from Equation 14 are smaller than the simulated timevalues. Thus, Equation 14 can be considered to be on the conservativeside. It is reliable to treat pressure data collected for t<t_(c) asthose corresponding to an infinitely long column. As a result, for a4-inch long shale sample, the no-flow boundary effect can be minimizedat X=1 inch. The valid time period for pressure measurement (0-3000 s)can cover a wide pressure range (100-750 psi) for the given example.

Numerical experiments are also conducted to check if the test proceduregives the “true” pressure-dependency of shale gas permeability. In anumerical experiment, the “true” permeability is that used as modelinput. Observed pressure data from the location about 1″ away from theinlet are used and random errors with magnitude of 0.2 psi are added tothe simulated pressures to consider the pressure measurement errors.FIG. 7 shows a graphic illustrating simulated gas pore gas pressure 700,real value 702 and with random error of +/−0.2 psi 704, as a function oftime at the core location about one inch away from the inlet, accordingto one example embodiment of the disclosure. The addition of randomerror does not make considerable difference in the pressure distributionbecause pressure measurement error is generally small.

As indicated by the line 602 in FIG. 6, the time when the no-flowboundary affects the pressure response in upstream locations increaseswith the distance to the no-flow boundary. So in order to ensure enoughtime for valid measurements, the pressure gauge should be put in areasonable distance away from the boundary. However at the same time, itshould not be too close to the inlet because the pressure response thereincreases from p_(i) to p₀ too fast. In the present method, themeasurement location is at X=1 inch.

While doing the actual measurements, an estimated time for validmeasurements can be calculated using Equation 14. The blue line 604(calculated from Equation 14) in FIG. 6 indicates that it is aconservative estimation and thus can be safely used in practice.

FIG. 8 shows a graphic illustrating comparison between the truepermeability as a function of pore gas pressure and the permeabilityfrom numerical experiment results, according to one example embodimentof the disclosure. As shown in the graph 800 in FIG. 8, results 804based on the laboratory test procedure discussed in the above sectionwith input k(p) and pressure data from numerical experiments are almostidentical to the “true” values 802 (or input k(p)), indicating that theproposed procedure is accurate and reliable. It can be observed thatthey are highly consistent with each other, which also means that therecorded pressure response is very close to that in the theoreticalmodel and the boundary effect is minimized at location X=1 inch.

Computer Readable Medium

In another example embodiment, the invention relates to computerprograms stored in computer readable media. Referring to FIG. 9, theforegoing process as explained with reference to FIGS. 1-8 can beembodied in computer-readable code. The code can be stored on, e.g., acomputer readable medium, such as a floppy disk 164, CD-ROM 162 or amagnetic (or other type) hard drive 160 forming part of a generalpurpose programmable computer. The computer, as known in the art,includes a central processing unit 150, a user input device such as akeyboard 154 and a user display 152 such as a flat panel LCD display orcathode ray tube display. According to this aspect of the invention, thecomputer readable medium includes logic operable to cause the computerto execute acts as set forth above and explained with respect to theprevious figures. The non-transitory computer-readable medium havingcomputer executable instructions cause a computer to perform theoperations of reading a measurement of a first pore gas pressure, p_(i),inside a sample assembly 300 comprising a sample of a subsurfaceformation, gas, and a pressure gauge. The instructions also includereading a measurement of a second pore gas pressure, p_(o), applied tothe inlet of a sample, where the second pore gas pressure is greaterthan the first pore gas pressure. The instructions also include readinga measurement of a third pore gas pressure, p, at location x at time tin the sample, and determining a total gas mass per unit volume of thesubsurface formation, m. The instructions also include determining apermeability of the subsurface formation, k, based at least in part onthe first pore gas pressure, the second pore gas pressure, the thirdpore gas pressure, and the gas density.

The computer executable instructions further cause the computer toperform the operation of determining the transport parameter of thesubsurface formation, D(p), using a first formula:

${D(p)} = \frac{\int_{p_{i}}^{p}{\frac{\lambda}{2}\frac{d\; m}{dp}{dp}}}{\frac{dp}{d\;\lambda}}$

where p_(i) is the first pore gas pressure inside the sample in assembly300 before the second pore gas pressure p_(o) is applied, p is the thirdpore gas pressure at location x at time t, m is the total gas mass perunit volume of the subsurface formation, and λ is an independentvariable calculated using the formula xt^(−1/2). Then permeability canbe determined from D(p) using Equation 9.

The computer executable instructions further cause the computer toperform the operation of determining the total gas mass per unit volumeof the subsurface formation, m, using a second formula:m=ϕρ+(1−ϕ)ρ_(a)

where ϕ is porosity of the subsurface formation, ρ is gas density of thenatural gas, and ρ_(a) is adsorbed gas mass per unit volume of thesubsurface formation.

The computer executable instructions further cause the computer toperform the operation of determining the porosity ϕ of the subsurfaceformation using a third formula:

$\phi = \frac{B - {A{\int_{p_{i}}^{p_{0}}{\lambda\frac{d\;\rho_{a}}{dp}{dp}}}}}{A{\int_{p_{i}}^{p_{0}}{\lambda\frac{d\;\left( {\rho - \rho_{a}} \right)}{dp}{dp}}}}$

where A is a cross-sectional area of the sample, and B is a slope of acurve of the cumulative gas flow into the sample at x=0 versus t^(1/2).

The computer executable instructions further cause the computer toperform the operation of determining the slope of the curve, B, using afourth formula:

$B = {A{\int_{p_{i}}^{p_{0}}{\lambda\frac{d\; m}{dp}{dp}}}}$

Methods according to the invention may provide improved estimates ofpermeability as a function of pore gas pressure and porosity ofsubsurface rock formations. Analytical models used to measurepressure-dependent gas permeability of shale are disclosed. Examplemethods and systems to measure shale gas permeability as a function ofpore gas pressure are disclosed. The advantages of new approach over thecurrently available ones include that it measures pressure-dependent gaspermeability more efficiently using a single test run and without anypresumption regarding a parametric relationship between gas permeabilityand pressure. In addition, the current invention also allows forestimating shale porosity from the related measurements.

While the invention has been described with respect to a limited numberof embodiments, those skilled in the art, having benefit of thisdisclosure, will appreciate that other embodiments can be devised whichdo not depart from the scope of the invention as disclosed herein.Accordingly, the scope of the invention should be limited only by theattached claims.

The invention claimed is:
 1. An apparatus for determining gaspermeability of a subsurface formation, comprising: a first sleevehaving a length, a diameter, a first open end, and a second open end; afirst end piece adapted to be inserted into the first open end and asecond end piece adapted to be inserted into the second open end of thefirst sleeve; a first hole formed through the first sleeve and into abody of a sample of the subsurface formation housed in the first sleeve;a half sleeve disposed on the first sleeve, the half sleeve having asecond hole corresponding to the first hole on the first sleeve; atubing comprising a pressure gauge, the tubing inserted through thefirst hole, the second hole, and into the body of the sample; ananchoring device for securing the tubing to the first sleeve and halfsleeve, thereby forming a sample assembly; and a pressure vessel forreceiving the sample assembly, the pressure vessel comprising at leasttwo ports connecting the sample to at least one pump with at least onepressure gauge, at least one port for applying a confining pressure, andat least one port for measuring pressure located at a known locationalong the sample, wherein determining gas permeability of the subsurfaceformation comprises: positioning the sample in the pressure vesselcomprising a fluid and the at least one pressure gauge; measuring afirst pore gas pressure, p₁, of a gas at a location, x, along thesample; applying a predetermined second pore gas pressure, p₀, to aninlet of the sample, the second pore gas pressure being greater than thefirst pore gas pressure; measuring a third pole gas pressure, p, as afunction of time, t, at a location, x, along the sample in the pressurevessel; in a computer, determining a gas density or total gas mass perunit volume of the subsurface formation, m; and in the computer,determining the gas permeability of the subsurface formation, k(p),based at least in part on the first pore gas pressure, the second poregas pressure, the third pore gas pressure as a function of time, and thegas density.
 2. The apparatus of claim 1, wherein the first sleeve andthe half sleeve comprise at least one of rubber and a polymericmaterial.
 3. The apparatus of claim 1, wherein an inner diameter of thehalf sleeve is smaller than the outer diameter of the first sleeve. 4.The apparatus of claim 1, wherein a length of the half sleeve is equalto or less than the length of the first sleeve.
 5. The apparatus ofclaim 1, wherein the fluid inside the pressure vessel comprises gas,water, or oil.
 6. The apparatus of claim 1, wherein the tubing furthercomprises a temperature gauge.
 7. The apparatus of claim 1, wherein thepressure gauge is coupled to the pressure tap using a flexible line. 8.The apparatus of claim 1, further comprising: an inlet pump configuredto pump gas from a first gas tank into the sample assembly inside thepressure vessel, an outlet pump configured to store gas from the sampleassembly inside the pressure vessel into a second gas tank, and gas inthe gas tanks, the sample, the pressure lines, and the pumps.
 9. Theapparatus of claim 1, wherein the step of determining a transportparameter of the subsurface formation further comprises: determining atransport parameter of the subsurface formation, D(p), using a firstformula:${D(p)} = \frac{\int_{p_{i}}^{p}{\frac{\lambda}{2}\frac{d\; m}{dp}{dp}}}{\frac{dp}{d\;\lambda}}$where λ is an independent variable calculated using the formulaλ=xt^(−1/2); and determining gas permeability k of the subsurfaceformation from D(p) using ${D(p)} = \frac{k\;\rho}{\mu}$ where μ standsfor gas viscosity, and ρ for gas density.
 10. The apparatus of claim 9,wherein the step of determining gas permeability of the subsurfaceformation further comprises: determining the total gas mass per unitvolume of the subsurface formation, m, using a second formula:m=ϕρ+(1−ϕ)ρ_(a) where ϕ is porosity of the subsurface formation, ρ isgas density of the gas, and ρ_(a) is adsorbed gas mass per unit volumeof the subsurface formation.
 11. The apparatus of claim 10, wherein thestep of determining gas permeability of the subsurface formation furthercomprises: determining the porosity 0 of the subsurface formation usinga third formula:$\phi = \frac{B - {A{\int_{p_{i}}^{p_{0}}{\lambda\frac{d\;\rho_{a}}{dp}{dp}}}}}{A{\int_{p_{i}}^{p_{0}}{\lambda\frac{d\;\left( {\rho - \rho_{a}} \right)}{dp}{dp}}}}$where A is a cross-sectional area of the sample, and B is a slope of acurve of the cumulative gas flow into the sample at x=0 versus t^(1/2).12. The apparatus of claim 11, wherein the step of determining gaspermeability of the subsurface formation further comprises: determiningthe slope of the curve, B, using a fourth formula:$B = {A{\int_{p_{i}}^{p_{0}}{\lambda\frac{d\; m}{dp}{{dp}.}}}}$
 13. Theapparatus of claim 1, wherein the subsurface formation comprises atleast one of shale, limestone, and sandstone.